\documentclass{article}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{amsfonts}
\usepackage{graphicx} 
\usepackage[left=1cm,right=1cm,top=1cm,bottom=2cm]{geometry}

\geometry{a4paper, margin=1in}

\title{Design and Implementation of Multigrid}
\author{3220101339 Jiang Zhou}
\date{2025/3/30}

\begin{document}

\maketitle


\section{Multigrid Solver Implementation}

\subsection{Class Overview}
The multigrid solver is implemented in C++ using the Eigen library for linear algebra operations. Two main classes are provided:

\begin{itemize}
    \item \texttt{Multigrid\_One\_Dimension} - Solves 1D Poisson-type equations
    \item \texttt{Multigrid\_Two\_Dimension} - Solves 2D Poisson-type equations
\end{itemize}

\subsection{Key Components}

\subsubsection{Data Structures}
\begin{itemize}
    \item \texttt{BoundaryConditions} struct - Handles different boundary condition types:
    \begin{itemize}
        \item Dirichlet
        \item Neumann
        \item Mixed boundary conditions
    \end{itemize}
    \item \texttt{MGParams} struct - Contains multigrid parameters:
    \begin{itemize}
        \item Grid resolution
        \item Restriction operator type (full weighting or injection)
        \item Interpolation type (linear or quadratic)
        \item Cycle type (V-cycle or FMG)
        \item Relaxation parameters
    \end{itemize}
\end{itemize}

\subsection{Grid Hierarchy}
The grid hierarchy is defined by:
\begin{equation}
h_\ell = \frac{1}{2^\ell}, \quad \ell = 0,1,\ldots,L
\end{equation}
where $L$ is the finest level.

\subsection{Smoothing Operation}
Weighted Jacobi iteration is used as the smoother:\(\omega = 0.75\)

\subsection{Intergrid Transfer Operators}
\subsubsection{Restriction (Fine-to-Coarse)}
\begin{itemize}
\item \textbf{Injection}: Direct transfer of values
\item \textbf{Full Weighting}: Weighted averaging of neighboring points
\end{itemize}

\subsection{Implementation Details}

\subsubsection{Transfer Operators}
The implementation provides flexible restriction and interpolation:

\begin{itemize}
    \item \textbf{Restriction}:
    \begin{itemize}
        \item Full weighting (averaging)
        \item Injection (simple copying)
    \end{itemize}
    
    \item \textbf{Interpolation}:
    \begin{itemize}
        \item Linear interpolation
        \item Quadratic interpolation
    \end{itemize}
\end{itemize}

\subsection{Definition: \(I_h^{2h}\) and \(I_{2h}^{h}\)}
The following C++ code implements the restriction and interpolation operations for a 2D multigrid method using the Eigen library. The code includes three main functions: \texttt{I\_h\_2h\_1d}, \texttt{I\_2h\_h\_1d}, and \texttt{I\_h\_2h} and \texttt{I\_2h\_h} for 2D operations. These functions handle different types of boundary conditions (Dirichlet, Neumann, and Mixed) and interpolation/restriction methods (linear and quadratic).

\subsubsection{1D Restriction: \texttt{I\_h\_2h\_1d}}
The function \texttt{I\_h\_2h\_1d} constructs a 1D restriction matrix from a fine grid to a coarse grid. The type of restriction is determined by the parameter \texttt{MG.restriction}.
\begin{itemize}
    \item \textbf{Quadratic Interpolation:} Transposes the restriction matrix.
    \item \textbf{Linear Interpolation:} Uses linear weights to interpolate coarse grid points to fine grid points.
\end{itemize}

\subsubsection{2D Restriction: \texttt{I\_h\_2h}}
The function \texttt{I\_h\_2h} constructs a 2D restriction matrix using the Kronecker product of 1D restriction matrices or custom logic for mixed boundary conditions.
\begin{itemize}
    \item \textbf{Injection:} Uses the Kronecker product of 1D injection matrices.
    \item \textbf{Full Weighting:} Uses the Kronecker product of 1D full weighting matrices or custom logic for mixed boundary conditions.
\end{itemize}

\subsubsection{2D Interpolation: \texttt{I\_2h\_h}}
The function \texttt{I\_2h\_h} constructs a 2D interpolation matrix using the Kronecker product of 1D interpolation matrices or custom logic for mixed boundary conditions.
\begin{itemize}
    \item \textbf{Quadratic Interpolation:} Transposes the 2D restriction matrix.
    \item \textbf{Linear Interpolation:} Uses custom logic to construct the interpolation matrix.
\end{itemize}

\subsection{\(I_h^{2h}\) and \(I_{2h}^{h}\)}
The following C++ code implements the restriction and interpolation operations for a 2D multigrid method using the Eigen library. The code includes three main functions: \texttt{I\_h\_2h\_1d}, \texttt{I\_2h\_h\_1d}, and \texttt{I\_h\_2h} and \texttt{I\_2h\_h} for 2D operations. These functions handle different types of boundary conditions (Dirichlet, Neumann, and Mixed) and interpolation/restriction methods (linear and quadratic).

\subsubsection{1D Restriction: \texttt{I\_h\_2h\_1d}}
The function \texttt{I\_h\_2h\_1d} constructs a 1D restriction matrix from a fine grid to a coarse grid. The type of restriction is determined by the parameter \texttt{MG.restriction}.

\begin{itemize}
    \item \textbf{Injection:} Directly maps fine grid points to coarse grid points. For example, if \(N = 8\), the injection matrix \(I_h^{2h}\) is:
    \[
    I_h^{2h} = 
    \begin{bmatrix}
    1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
    0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\
    0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\
    0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \\
    \end{bmatrix}
    \]
    \item \textbf{Full Weighting:} Uses weights to average neighboring fine grid points to form the coarse grid points. For example, if \(N = 8\), the full weighting matrix \(I_h^{2h}\) is:
    \[
    I_h^{2h} = 
    \begin{bmatrix}
    1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
    0.25 & 0.5 & 0.25 & 0 & 0 & 0 & 0 & 0 & 0 \\
    0 & 0 & 0 & 0.25 & 0.5 & 0.25 & 0 & 0 & 0 \\
    0 & 0 & 0 & 0 & 0 & 0 & 0.25 & 0.5 & 0.25 \\
    \end{bmatrix}
    \]
\end{itemize}

\subsubsection{2D Restriction: \texttt{I\_h\_2h}}
The function \texttt{I\_h\_2h} constructs a 2D restriction matrix using the Kronecker product of 1D restriction matrices or custom logic for mixed boundary conditions.

\begin{itemize}
    \item \textbf{Injection:} Uses the Kronecker product of 1D injection matrices. The 2D injection matrix \(I_h^{2h}\) is:
    \[
    I_h^{2h} = I_h^{2h}(1d) \otimes I_h^{2h}(1d) 
    \]
    \item \textbf{Full Weighting:} Uses the Kronecker product of 1D full weighting matrices or custom logic for mixed boundary conditions. For example, if \(N = 8\), the 2D full weighting matrix \(I_h^{2h}\) is:
    \[
    I_h^{2h} = I_h^{2h}(1d) \otimes I_h^{2h}(1d) =
    \begin{bmatrix}
    1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
    0 & 0.25 & 0.5 & 0.25 & 0 & 0 & 0 & 0 & 0 \\
    0 & 0 & 0 & 0.25 & 0.5 & 0.25 & 0 & 0 & 0 \\
    0 & 0 & 0 & 0 & 0 & 0.25 & 0.5 & 0.25 & 0 \\
    0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1\\
    \end{bmatrix}
    \otimes
    I_h^{2h}(1d)
    \]
\end{itemize}

\subsubsection{2D Interpolation: \texttt{I\_2h\_h}}
The function \texttt{I\_2h\_h} constructs a 2D interpolation matrix using the Kronecker product of 1D interpolation matrices or custom logic for mixed boundary conditions.

\begin{itemize}
    \item \textbf{Quadratic Interpolation:} Transposes the 2D restriction matrix. For example, if \(N = 8\), the 2D quadratic interpolation matrix \(I_{2h}^{h}\) is:
    \[
    I_{2h}^{h} = (I_h^{2h}(1d) \otimes I_h^{2h}(1d))^T
    \]
    \item \textbf{Linear Interpolation:} Uses custom logic to construct the interpolation matrix.
\end{itemize}

\subsection{Boundary Conditions}
Implemented types:
\begin{itemize}
\item Dirichlet: $u|_{\partial\Omega} = g$
\item Neumann: $\frac{\partial u}{\partial n}|_{\partial\Omega} = h$
\item Mixed: Combination of both
\end{itemize}

\subsection{Discretization}
The Laplacian is discretized using finite differences:
\begin{equation}
\Delta u \approx \frac{u_{i-1} - 2u_i + u_{i+1}}{h^2} \quad \text{(1D)}
\end{equation}
\begin{equation}
\Delta u \approx \frac{u_{i-1,j} + u_{i+1,j} + u_{i,j-1} + u_{i,j+1} - 4u_{i,j}}{h^2} \quad \text{(2D)}
\end{equation}

\subsection{Performance Features}

\begin{itemize}
    \item Sparse matrix storage using Eigen's \texttt{SparseMatrix}
    \item Recursive grid coarsening
    \item Convergence monitoring with relative residual checks
    \item Timing measurements for performance analysis
    \item Error computation and convergence rate tracking
\end{itemize}

\subsection{Usage Example}

The solver is configured through JSON input files specifying.\\
If the input.json file detects an invalid input, the terminal will report an error and give a reasonable input. \\
If I give an undefined funtion like \(\sin(x)*\cos(y)\), then you can get :
\begin{figure}
    \centering
        \includegraphics[width=0.5\linewidth]{figure/input.png}
        \caption{Invaild Input}
        \label{input}
\end{figure}

\section{Problem 2}
\subsection{\(u(x,y) = e^{y+\sin x}\)}
For Function1(\(u(x,y) = e^{y+\sin x}\)) and Function1-1d(\(u(x) = e^{x+\sin x}\)), we set \(\epsilon = 10^{-8}\) and the zero-vector initial guess:

\begin{figure}
    \centering
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_-(-sin(x) + cos^2(x) + 1) * exp(y + sin(x))_Mixed_1e-08_full_weighting_quadratic.png}
        \caption{Function1:Mixed V-Cycles full weighting quadratic}
        \label{function1_1}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_-(-sin(x) + cos^2(x) + 1) * exp(y + sin(x))_NEUMANN_1e-08_full_weighting_quadratic.png}
        \caption{Function1:NEUMANN V-Cycles full weighting quadratic}
        \label{function1_2}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_-(-sin(x) + cos^2(x) + 1) * exp(y + sin(x))_DIRICHLET_1e-08_full_weighting_quadratic.png}
        \caption{Function1:Dirichlet V-Cycles full weighting quadratic}
        \label{function1_3}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_-(-sin(x) + cos^2(x) + 1) * exp(y + sin(x))_DIRICHLET_1e-08_injection_linear.png}
        \caption{Function1:Dirichlet V-Cycles injection linear}
        \label{function1_4}
    \end{minipage}\hfill
\end{figure}

\begin{figure}
    \centering
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_-(-sin(x) + cos^2(x) + 1) * exp(y + sin(x))_NEUMANN_1e-08_injection_linear.png}
        \caption{Function1:Mixed V-Cycles injection linear}
        \label{function1_5}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-sin(x) + cos^2(x) + 1) * exp(x + sin(x))_V_Cycles_DIRICHLET_1e-08_full_weighting_quadratic.png}
        \caption{Function1 1d:DIRICHLET V-Cycles full weighting quadratic}
        \label{function1_6}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-sin(x) + cos^2(x) + 1) * exp(x + sin(x))_V_Cycles_Mixed_1e-08_full_weighting_quadratic.png}
        \caption{Function1 1d:Mixed V-Cycles full weighting quadratic}
        \label{function1_7}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-sin(x) + cos^2(x) + 1) * exp(x + sin(x))_FMG_Mixed_1e-08_full_weighting_quadratic.png}
        \caption{Function1 1d:Mixed FMG full weighting quadratic}
        \label{function1_8}
    \end{minipage}\hfill
\end{figure}

\begin{figure}
    \centering
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-sin(x) + cos^2(x) + 1) * exp(x + sin(x))_V_Cycles_DIRICHLET_1e-08_injection_linear.png}
        \caption{Function1 1d:DIRICHLET V-Cycles injection linear}
        \label{function1_9}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-sin(x) + cos^2(x) + 1) * exp(x + sin(x))_FMG_DIRICHLET_1e-08_full_weighting_quadratic.png}
        \caption{Function1 1d:DIRICHLET FMG full weighting quadratic}
        \label{function1_10}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-sin(x) + cos^2(x) + 1) * exp(x + sin(x))_V_Cycles_DIRICHLET_1e-08_full_weighting_linear.png}
        \caption{Function1 1d:DIRICHLET V-Cycles full weighting quadratic}
        \label{function1_11}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-sin(x) + cos^2(x) + 1) * exp(x + sin(x))_FMG_DIRICHLET_1e-08_full_weighting_linear.png}
        \caption{Function1 1d:DIRICHLET FMG full weighting quadratic}
        \label{function1_12}
    \end{minipage}\hfill
\end{figure}

\subsection{\(u(x,y) =\sin(\pi x) \cdot \sin(\pi y)\)}
For Function1(\(u(x,y) =\sin(\pi x) \cdot \sin(\pi y)\)) and Function1-1d(\(u(x) =\sin^2(\pi x)\)), we set \(\epsilon = 10^{-8}\) and the zero-vector initial guess:
\begin{figure}
    \centering
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_DIRICHLET_1e-08_full_weighting_quadratic.png}
        \caption{Function2:DIRICHLET V-Cycles full weighting quadratic}
        \label{function2_1}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_DIRICHLET_1e-08_full_weighting_linear.png}
        \caption{Function2:DIRICHLET V-Cycles full weighting linear}
        \label{function2_2}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_FMG_DIRICHLET_1e-08_full_weighting_linear.png}
        \caption{Function2:DIRICHLET FMG full weighting linear}
        \label{function2_3}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_FMG_DIRICHLET_1e-08_injection_linear.png}
        \caption{Function2:DIRICHLET FMG injection linear}
        \label{function2_4}
    \end{minipage}\hfill
\end{figure}

\begin{figure}
    \centering
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_V_Cycles_NEUMANN_1e-08_injection_linear.png}
        \caption{Function2:NEUMANN V-Cycles injection linear}
        \label{function2_5}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_V_Cycles_NEUMANN_1e-08_full_weighting_linear.png}
        \caption{Function2:NEUMANN V-Cycles full weighting linear}
        \label{function2_6}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_FMG_NEUMANN_1e-08_full_weighting_linear.png}
        \caption{Function2:NEUMANN FMG full weighting linear}
        \label{function2_7}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_FMG_NEUMANN_1e-08_full_weighting_quadratic.png}
        \caption{Function2:NEUMANN FMG full weighting quadratic}
        \label{function2_8}
    \end{minipage}\hfill
\end{figure}

\begin{figure}
    \centering
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_V_Cycles_NEUMANN_1e-08_full_weighting_linear.png}
        \caption{Function2:NEUMANN V-Cycles full weighting linear}
        \label{function2_9}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_V_Cycles_Mixed_1e-08_full_weighting_quadratic.png}
        \caption{Function2:Mixed V-Cycles full weighting quadratic}
        \label{function2_10}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_FMG_Mixed_1e-08_full_weighting_quadratic.png}
        \caption{Function2:Mixed FMG full weighting quadratic}
        \label{function2_11}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_2 * PI^2 * sin(PI * x) * sin(PI * y)_V_Cycles_Mixed_1e-08_injection_linear.png}
        \caption{Function2:Mixed V-Cycles injection linear}
        \label{function2_12}
    \end{minipage}\hfill
\end{figure}

\begin{figure}
    \centering
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_2 * PI^2 * sin(PI * x) * sin(PI * x)_V_Cycles_DIRICHLET_1e-08_full_weighting_linear.png}
        \caption{Function2-2d:DIRICHLET V-Cycles full weighting linear}
        \label{function2_13}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_2 * PI^2 * sin(PI * x) * sin(PI * x)_V_Cycles_DIRICHLET_1e-08_full_weighting_quadratic.png}
        \caption{Function2-2d:DIRICHLET V-Cycles full weighting quadratic}
        \label{function2_14}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_2 * PI^2 * sin(PI * x) * sin(PI * x)_V_Cycles_NEUMANN_1e-08_full_weighting_quadratic.png}
        \caption{Function2-2d:NEUMANN FMG full weighting quadratic}
        \label{function2_15}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_2 * PI^2 * sin(PI * x) * sin(PI * x)_FMG_DIRICHLET_1e-08_full_weighting_quadratic.png}
        \caption{Function2-2d:DIRICHLET FMG full weighting quadratic}
        \label{function2_16}
    \end{minipage}\hfill
\end{figure}

\subsection{\(u(x,y) = e^{\sin(x) + \cos(y)}\)}
For Function3(\(u(x,y) = e^{\sin(x) + \cos(y)}\)) and Function3-1d(\(u(x) = e^{\sin(x) + \cos(x)}\)), we set \(\epsilon = 10^{-8}\) and the zero-vector initial guess:

\begin{figure}
    \centering
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_-(-cos(y) - sin(x) + sin^2(y) + cos^2(x)) * exp(sin(x) + cos(y))_V_Cycles_DIRICHLET_1e-08_full_weighting_linear.png}
        \caption{Function3:DIRICHLET V-Cycles full weighting linear}
        \label{function3_1}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_-(-cos(y) - sin(x) + sin^2(y) + cos^2(x)) * exp(sin(x) + cos(y))_V_Cycles_DIRICHLET_1e-08_full_weighting_quadratic.png}
        \caption{Function3:DIRICHLET V-Cycles full weighting quadratic}
        \label{function3_2}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_-(-cos(y) - sin(x) + sin^2(y) + cos^2(x)) * exp(sin(x) + cos(y))_FMG_DIRICHLET_1e-08_full_weighting_quadratic.png}
        \caption{Function3:DIRICHLET FMG full weighting quadratic}
        \label{function3_3}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_2d_-(-cos(y) - sin(x) + sin^2(y) + cos^2(x)) * exp(sin(x) + cos(y))_V_Cycles_DIRICHLET_1e-08_injection_linear.png}
        \caption{Function3:DIRICHLET V-Cycles injection linear}
        \label{function3_4}
    \end{minipage}\hfill
\end{figure}


\begin{figure}
    \centering
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-cos(x) - sin(x) + sin^2(x) + cos^2(x)) * exp(sin(x) + cos(x))_V_Cycles_DIRICHLET_1e-08_injection_linear.png}
        \caption{Function3-1d:DIRICHLET V-Cycles injection linear}
        \label{function3_5}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-cos(x) - sin(x) + sin^2(x) + cos^2(x)) * exp(sin(x) + cos(x))_FMG_DIRICHLET_1e-08_injection_linear.png}
        \caption{Function3-1d:DIRICHLET FMG injection linear}
        \label{function3_6}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-cos(x) - sin(x) + sin^2(x) + cos^2(x)) * exp(sin(x) + cos(x))_FMG_DIRICHLET_1e-08_full_weighting_linear.png}
        \caption{Function3-1d:DIRICHLET FMG full weighting linear}
        \label{function3_7}
    \end{minipage}\hfill
    \begin{minipage}{0.23\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/test_1d_-(-cos(x) - sin(x) + sin^2(x) + cos^2(x)) * exp(sin(x) + cos(x))_FMG_DIRICHLET_1e-08_full_weighting_linear.png}
        \caption{Function3-1d:DIRICHLET V-Cycles full weighting linear}
        \label{function3_8}
    \end{minipage}\hfill
\end{figure}

\section{Problem 3}
To determine the critical value of $\epsilon$ where the multigrid solver fails to achieve the preset accuracy, we conducted a series of experiments with gradually decreasing $\epsilon$ values. The results are summarized below:
The program begins to fail consistently when $\epsilon$ reaches approximately \boxed{2.3 \times 10^{-16}}. \\
We change the MultigridTwoDimension.solve() in order to solve the question, as shown in the picture:
\begin{figure}
    \centering
    \includegraphics[width=0.5\linewidth]{figure/problem3.png}
    \caption{Problem3: code}
    \label{Problem3}
\end{figure}
We test on the function2 with the initialization below:
\begin{itemize}
    \item grid resolution: $8$
    \item restriction full: weighting
    \item interpolation: quadratic
    \item BoundaryConditions: DIRICHLET
\end{itemize}
Then, we can get the reason of the critical value of $\epsilon$:
\begin{itemize}
\item The residual norm stagnates around machine epsilon ($\approx 2.3 \times 10^{-16}$)
\item Further iterations don't improve the solution quality
\end{itemize}

\section{B}
In order to compare multigrid with Poisson Solver, we can set the "recursion limit": 0 and "max iterations": 1 in order to use Poisson Solver to solve the problem.\\
Then we can get the result below:
\begin{figure}
    \centering
    \begin{minipage}{0.45\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/B1.png}
        \caption{B: Multigrid}
        \label{B1}
    \end{minipage}\hfill
    \begin{minipage}{0.45\linewidth}
        \centering
        \includegraphics[width=\linewidth]{figure/B2.png}
        \caption{B: Poisson Solver}
        \label{B2}
    \end{minipage}\hfill
\end{figure}
\end{document}